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Abstract 

A model accounting for finite spatial dimensions of the deposit patterns in the evap- 
orating sessile drops of colloidal solution on a plane substrate is proposed. The model is 
based on the assumption that the solute particles occupy finite volume and hence these 
dimensions are of the steric origin. Within this model, the geometrical characteristics 
of the deposition patterns are found as functions of the initial concentration of the 
solute, the initial geometry of the drop, and the time elapsed from the beginning of the 
drying process. The model is solved analytically for small initial concentrations of the 
solute and numerically for arbitrary initial concentrations of the solute. The agreement 
between our theoretical results and the experimental data is demonstrated, and it is 
shown that the observed dependence of the deposit dimensions on the experimental pa- 
rameters can indeed be attributed to the finite dimensions of the solute particles. These 
results are universal and do not depend on any free or fitting parameters; they are im- 
portant for understanding the evaporative deposition and may be useful for creating 
controlled deposition patterns. 

PACS: 47.55. Dz — Drops and bubbles; 68.03.Fg — Evaporation and condensation; 
81.15.-z — Methods of deposition of films and coatings; film growth and epitaxy. 

1 Introduction 

The problem of the so-called "coffee-drop deposit" has recently aroused great interest. The 
residue left when coffee dries on the countertop is usually darkest and hence most concentrated 
along the perimeter of the stain. Ring-like stains, with the solute segregated to the edge of 
a drying drop, are not particular to coffee. Mineral rings left on washed glassware, banded 
deposits of salt on the sidewalk during winter, and enhanced edges in water color paintings 
are all examples of the variety of physical systems displaying similar behavior and understood 
by coffee-drop deposit terminology. 



* Current address: Department of Physics, University of Michigan, 500 E. University Ave., Ann Arbor, 
MI 48109. E-mail: yopopov@umich.edu 



1 



Understanding the process of drying of such solutions is important for many scientific 
and industrial applications, where ability to control the distribution of the solute during 
drying process is at stake. For instance, in the paint industry, the pigment should be evenly 
dispersed after drying, and the segregation effects are highly undesirable. Also, in the pro- 
tein crystallography, attempts are made to assemble the two-dimensional crystals by using 
evaporation driven convection PflElEI, and hence solute concentration gradients should be 
avoided. On the other hand, in the production of nanowires ^ or in surface patterning [5] 
perimeter-concentrated deposits may be of advantage. Recent important applications of this 
phenomenon related to DNA stretching in a flow have emerged as well [6J. For instance, 
a high-throughput automatic DNA mapping was suggested [Zj, where fluid flow induced by 
evaporation is used for both stretching DNA molecules and depositing them onto a sub- 
strate. Droplet drying is also important in the attempts to create arrays of DNA spots for 
gene expression analysis. 

Ring-like deposit patterns have been studied experimentally by a number of groups. Dif- 
ficulties of obtaining a uniform deposit jH], deformation of sessile drops due to a sol-gel 
transition of the solute at the contact line [HI E]], stick-slip motion of the contact line of 
colloidal liquids ^], multiple ring formation [H], and the effect of ring formation on 
the evaporation of the sessile drops were all reported. The evaporation of the sessile 
drops (regardless of the solute presence) has also been investigated extensively. Constancy 
of the evaporation flux was demonstrated ^3 ^| , and the change of the geometrical char- 
acteristics (contact angle, drop height, contact-line radius) during drying was measured in 
detail |I2l [TBI CHI ED! ■ 

The most recent and complete experimental effort to date on coffee-drop deposits was 
conducted by Robert Deegan et al. j^H 1221 I2H1 121] ■ Most experimental data referred to 
in this work originate from observations and measurements of this group. They reported 
extensive results on ring formation and demonstrated that these could be quantitatively 
accounted for. The main ideas of the theory of solute transfer in such physical systems have 
also been developed in their work |2I]. It was observed that the contact line of a drop of 
liquid remains pinned during most of the drying process. While the highest evaporation 
occurs at the edges, the bulk of the solvent is concentrated closer to the center of the drop. 
In order to replenish the liquid removed by evaporation at the edge, a flow from the inner to 
the outer regions must exist inside the drop. This flow is capable of transferring all of the 
solute to the contact line and thus accounts for the strong contact-line concentration of the 
residue left after complete drying. This theory is very robust since it is independent of the 
nature of the solute and only requires pinning of the edge during drying (which can occur in 
a number of possible ways: surface roughness, chemical heterogeneities etc). Among other 
things, we will reproduce some of its results in this work. 

Mathematically, the most complicated task is related to determining the evaporation rate 
from the surface of the drop. An analogy between the diffusive concentration fields and 
the electrostatic potential fields was suggested |23 12^] , so that an equivalent electrostatic 
problem can be solved instead of the evaporation problem. Important analytical solutions to 
this equivalent problem in various geometries were first derived by Lebedev |2^, and a few 
useful consequences from these analytical results were later reported in Ref. 

In this work, we continue development of the theory of solute transfer and deposit growth. 
Most previous works address the issue of the deposit mass accumulation at the drop boundary. 
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however, they treat the solute particles as if they do not occupy any volume, and hence all 
the solute can be accommodated at the one-dimensional singularity of the contact line. In 
reality, the solute deposit accumulated at the perimeter has some thickness, and the shape 
of the residue in a round drop resembles a ring rather than an infinitely thin circumference 
of the circle. The earlier efforts were aimed at describing how the mass of the contact-line 
deposit grows with time and how it depends on such geometrical characteristics of the drop 
as its radius (for circular drops [211 122]) or its opening angle and the distance from the 
vertex (for pointed drops |28,,2n|)- Little attempt has been made to describe the geometrical 
characteristics of the contact-line deposit itself, for instance, the width and the height of the 
deposit ring. At the same time, there is solid experimental data [2311211 on various geometrical 
characteristics of the ring and their dependence on time, the initial solute concentration, and 
the drop geometry. Here we develop a simple model that addresses this lack of understanding 
of the geometrical properties of the contact-line deposit and accounts for the finite size of 
the deposit ring. We attribute the finite volume of the deposit simply to the finite size of 
the solute particles, i.e. we assume the particles do occupy some volume and hence cannot 
be packed denser than certain concentration. The model is solved in the simplest case of the 
circular geometry both analytically and numerically, and the results of the two methods are 
compared with the experimental data of Refs. [2S1 121] (and with each other). It turns out 
that this model is sufficient to explain most of the collected data. It should be noted that 
the model is as universal and robust (in its range of validity) as the zero-volume theory of 
Deegan et al. [211 122] since it is based on essentially the same physical principles. 

The notion that the profile of the deposit could be found by the simple assumption that 
the solute becomes immobilized when the volume fraction reaches a threshold was originally 
suggested by Todd Dupont [SU]. First efforts to create a model were conducted by Robert 
Deegan [^ |21] who formulated some physical assumptions, wrote them down mathemati- 
cally, and obtained some early-time exponents. Here we present the entire problem, including 
its full formulation and its analytical and numerical solutions (not reported previously). In 
the next section, we formulate the model, describe the system, and address some issues of 
the geometry and the evaporation rate. Then, we derive the governing equations from the 
conservation of mass and later solve them analytically for small initial concentrations of the 
solute and numerically for arbitrary initial concentrations of the solute. A discussion section 
concludes this work. 

2 Model, assumptions, and geometry 

System. We consider a sessile droplet of solution on a horizontal surface (substrate). 
The nature of the solute is not essential for the mechanism. The typical diameter of the 
solute particles in Deegan's experiments [211 1221I2S1I21] was of the order of 0.1-1 yum; we will 
assume a similar order of magnitude throughout this work. For smaller particles diffusion 
becomes important; for larger particles sedimentation may play an important role. 

The droplet is bounded by the contact line in the plane of the substrate. This (macro- 
scopic) contact line is defined as the common one- dimensional boundary of all three phases 
(liquid, air and solid substrate). We will restrict our attention to the case of the round drops, 
which is both of most practical importance and the easiest to treat mathematically. 
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We assume that the droplet is sufficiently small so that the surface tension is dominant, 
and the gravitational effects can be neglected. Mathematically, the balance of the gravita- 
tional force and the surface tension is controlled by the ratio of the (maximal) hydrostatic 
pressure pghmax to the Laplace pressure 2ahmax/R'i, where p is the fluid density, g is the 
gravitational constant, a is the surface tension at the liquid- air interface, Ri is the drop radius 
in the plane of the substrate, and hmax is the maximal height of the drop. For the typical 
experimental conditions this ratio pgRf/2a is quite small (about 0.25), and thus gravity is 
indeed unimportant and the surface shape is governed mostly by the surface tension. Our 
treatment will produce the main-order term in the expansion in this parameter, and since 
the parameter value is not an order of magnitude smaller than one, it may be necessary to 
construct the correctional terms for better quantitative agreement. For the present purposes, 
even the main term turns out to be sufficient to obtain the agreement with the experimental 
results. 

Experimentally, the contact line remains pinned during most of the drying process. There- 
fore, we do not assume that the contact angle 6 between the liquid-air interface and the plane 
of the substrate is constant in time. A strongly pinned contact line can sustain a wide range 
of (macroscopic) contact angles. The pinning mechanism can be described as self-pinning, 
i.e. pinning by the deposit brought to the contact line by the hydrodynamic flows caused by 
evaporation. A pinned contact line entails fluid flow toward that contact line. The "elasticity" 
of the liquid-air interface fixed at the contact line provides the force driving this flow. 

We will deal with small contact angles (6' ^ 1) as it is almost always the case in the 
experimental reahzations, including the experiments of Ref. ^^2^ (typically, 9max < 0.1-0.3). 
It will also be seen necessary to assume that the contact angle is small in order to obtain any 
analytical results in a closed form. A drop with a small contact angle is necessarily thin, i.e. 
its maximal height is much smaller than its radius and the slope of the free surface is small 
(|V/i| <C 1). Thus, we consider small contact angles, or, equivalently, thin drops. 

We also consider slow flows, i.e. flows with low Reynolds numbers, which amounts to the 
neglect of the inertial terms in the Navier-Stokes equation. 

The free surface is described by the local mean curvature that is spatially uniform at any 
given moment of time, but changes with time as the droplet dries. Ideally, the surface shape 
should be considered dynamically together with the flow field inside the drop. However, 
as was shown earlier jSHl 1^ , for fiow velocities much lower than the characteristic velocity 
V* = a/Sr] (where a is the surface tension and r] is the dynamic viscosity), which is about 
24 m/s for water under normal conditions, one can consider the surface shape independently 
of the fiow and use the equilibrium result at any given moment of time for finding the fiow at 
that time. Equivalently, the ratio of the viscous forces to the capillary forces is the capillary 
number Ca = rjv/a (where v is some characteristic value of the fiow velocity, which is of the 
order of 1-10 /im/s), and this number is of the order of 10~^-10~^ under typical experimental 
conditions. Thus, the capillary forces are by far the dominant ones. 

Geometry and surface shape. The cylindrical coordinates (r, (j), z) will be used through- 
out this work, as they are most natural for the geometry of interest. The origin is chosen 
in the center of the circular footprint of the drop on the substrate. Coordinate z is always 
normal to the substrate, and the substrate is described by z = 0, with z being positive on the 
droplet side of the space. Coordinates (r, 0) are the polar radius and the azimuthal angle. 
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respectively, so that the contact hne is described by r = where Ri is the radius of the 
drop footprint. Due to the axial symmetry of the problem and our choice of the coordinates, 
no quantity depends on the azimuthal angle 0. 

Our model pictures the drop as a two-component system (the components being "the 
fluid" and "the solute"), which has two "phases": "the liquid phase" in the middle of the 
drop and "the deposit phase" near the contact line. Both components are present in both 
phases, and the difference between the phases lies only in the concentration of the solute in 
each phase. In the deposit phase, the volume fraction of the solute p is high and fixed in both 
space and time. Thus, p is just a constant number, one can think of it as comparable to the 
close-packing fraction or unity. (The case of p = 1 may seem to be special as there is no fluid 
in the deposit phase; however, for small initial concentrations of the solute this case will be 
seen to lead to exactly the same main order results.) In the liquid phase, the volume fraction 
of the solute x varies in space and changes with time, and it is relatively small compared to 
p. The initial volume fraction Xi = x(0) is constant throughout the drop; at later moments 
the solute gets redistributed due to the flows, and the concentration becomes different in 
different parts of the liquid phase. The volume fraction of the fluid is then {1 — p) in the 
deposit phase and (1 — x) in the liquid phase. Note that we do not require Xj <^ p so far, 
although we do assume Xi < P- It should also be emphasized that we do not presume there is 
any real "phase difference" between the so-called phases: one phase is just defined as having 
the maximal reachable solute fraction p (the solute cannot move in this phase) while the 
other phase is characterized by lower solute fraction x (in this phase the solute can move and 
hence its concentration can change in time and space). Besides this difference, the phases 
are essentially identical. The idea that the solute loses its mobility when its concentration 
exceeds some threshold was suggested by Todd Dupont pUj . 

Since the drop is thin, we employ the vertically averaged flow velocity 



where u<.(r, z) is the in-plane component of the local three-dimensional velocity u{r,z), and 
ht{r) is the thickness of the drop at distance r from the center. By making this approximation, 
we implicitly assume that there is no vertical segregation of the solute, and thus we turn our 
model into effectively two-dimensional. This is done mostly for simplicity and is not expected 
to affect our main conclusions (see Discussion). Within this model, it is natural to assume 
that the boundary between the phases is vertical. Thus, the particles get stacked uniformly 
at all heights when being brought to the phase boundary by the hydrodynamic flow v(r). 
This boundary can be pictured as a vertical wall at some radius R{t) from the center of the 
drop, and this wall propagates from the contact line [located at Ri = -R(O)] towards the center 
of the drop. Fig. ^illustrates the mutual location of the two phases, and Fig. |21 schematically 
shows the time evolution of the drying process and growth of the deposit phase. 

The geometrical parameters of the model are shown in Fig. |3| The radius of the drop 
is Ri, the radius of the phase boundary is R{t), and -R(O) = Ri. The height of the phase 
boundary is H{t), and the initial condition is H{0) = 0. In the liquid phase, we conveniently 
split the total height of the free surface ht{r,t) into the sum of H{t) and h{r,t). 

Since H is independent of r, function h{r) satisfies the Young-Laplace equation (the 




(1) 
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Figure 1: Mutual location of the two "phases" in the drying drop: L is the liquid phase, and 
D is the deposit phase. 




Figure 2: Time evolution of the deposit phase growth: side view (left) and top view (right). 
Only the deposit phase is shown. Thickness of the ring is exaggerated compared to the typical 
experimental results. 




Hr,t) 




R(t) 



Figure 3: Geometry of the problem. Vertical scale is exaggerated in order to see the details, 
typically H <ti Ri and h <ti Ri. 
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statement of the mechanical equihbrium of the hquid-air interface) 

2K = (2) 

(7 

where a is the surface tension, Ap is the pressure difference across the hquid-air interface, 
and K is the mean curvature of the surface, uniquely related to the surface shape h by 
differential geometry. For typical drying conditions, Ap and h vary with time slowly. As 
was shown earlier [211 EI], it is sufficient to find the equilibrium surface shape first, and 
then determine the velocity field for this fixed functional form of h with time being just an 
adiabatic parameter, instead of solving for all the dynamical quantities simultaneously. Thus, 
the right-hand side of Eq. Q does not depend on the local coordinates of a point within the 
drop (although it does depend on time), and the equation itself expresses the global condition 
of spatial constancy of the mean curvature throughout the drop. It defines the equilibrium 
surface shape at any given moment of time. 

The solution to Eq. (0) with boundary condition h{R) = is just a spherical cap, and 
hence the shape of the upper part of the drop (above the dashed line in Fig. E)) is just a 
spherical cap: 



h{r,t) 



^'^^^ r2-i?(t) cot 0(t). (3) 



\ sin2^(t) 



Here 9{t) is the angle between the liquid-air interface and the substrate at phase boundary, 
and functions R{t) and 0{t) are related via the right-hand side of Eq. 

^W = ^«in^(t). (4) 

In the limit of small contact angles, 6' <^ 1, the preceding expression adopts even simpler 
form: 



r 



2 



KT,t) = ^-^^^^e{t) + o{e^). (5) 

Note that we do not assume that 9{t) and h{r, t) are necessarily positive at all times: both can 
be negative at later drying stages, and the shape of the liquid-air interface may be concave. 
Both convex and concave solutions for h{r,t) are consistent with Eq. the right-hand side 
of this equation can have either sign. By definition, both 6{t) and h{r,t) are positive when 
the surface is convex (and hence they are positive at the beginning of the drying process) 
and negative when the surface is concave. The initial value of 6{t) coincides with the initial 
contact angle 9i = 9{0). 

Clearly, there are three unknown functions of time in this geometry: 6{t), R{t) and H{t). 
However, these quantities are not independent of each other. Since we assume that the solute 
particles fill up the entire space between the substrate and the liquid-air interface when being 
brought to the phase boundary, the three geometrical functions are related by the constraint 

Physically, the angles between the liquid-air interface and the substrate are identical on 
both sides of the phase boundary [6 = \dH/dR\), and hence h{r) and its first derivative are 



7 



continuous past this boundary. Thus, there are actually only two independent functions of 
time, 9{t) and R{t). Condition (0) was first introduced by Robert Deegan 1^. 

The geometrical definitions above allow one to determine the volume of each of the two 
phases. Volume of the liquid phase is simply 

Vl = I ihir,t) + Hit)) 2nrdr = 2n + + 0(6'). (7) 

Taking into account relation (jU)), an infinitesimal variation of this volume can be expressed 
via the infinitesimal variations of 6 and R: 

dVL = ^d + 2'KHRdR. (8) 

The first term is responsible for the motion of the liquid-air interface, and the second term 
corresponds to the inward shift of the phase boundary. It is also straightforward to obtain 
an expression for the differential of the volume of the deposit phase, which has only the term 
related to the inward shift of the phase boundary: 

dVD = -27iHRdR. (9) 

We will use the last two expressions in the following section. We will also adopt the notation 
that subscripts L and D refer to the liquid and deposit phases, respectively. 

Before proceeding to the main section, we will make a note on the evaporation rate. 



Evaporation rate. In order to determine the flow caused by evaporation, one needs to 
know the flux profile of liquid leaving each point of the surface. This quantity is independent 
of the processes going on inside the drop and must be determined prior to considering any 
such processes. 

The functional form of the evaporation rate J(r) (defined as the evaporative mass loss per 
unit surface area per unit time) depends on the rate-limiting step, which can, in principle, 
be either the transfer rate across the liquid- vapor interface or the diffusive relaxation of the 
saturated vapor layer immediately above the drop. We assume that the rate-limiting step 
is the diffusion of the saturated vapor. Indeed, the transfer rate across the liquid-vapor 
interface is characterized by the time scale of the order of 10~^° s, while the diffusion process 
has characteristic times of the order of R^/D (where D is the diffusion constant for vapor 
in air and Ri is a characteristic size of the drop), which is of the order of seconds for water 
drops under typical drying conditions. The diffusion-limited evaporation rapidly attains a 
steady state. Indeed, the ratio of the time required for the vapor-phase water concentration 
to adjust to the changes in the droplet shape {R? /D) to the droplet evaporation time is of 
the order of (n^ — rioo)! p ~ 10~^, where is the density of the saturated vapor just above 
the liquid-air interface, rzoo is the ambient vapor density, and p is the fluid density PZj, i.e. 
the vapor concentration adjusts rapidly compared to the evaporation time. 

As the rate-limiting process is the diffusion, vapor density n above the liquid-vapor inter- 
face obeys the diffusion equation. Since the process is quasi-steady, this diffusion equation 
reduces to the Laplace equation 

V^n = 0. (10) 
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in substrate 

Figure 4: Illustration of the analogy between the evaporation rate J for a liquid drop and the 
electric field E for a conductor. Consideration of the drop (or conductor) and its reflection 
in the plane of the substrate significantly simplifies the boundary problem. 



This equation is to be solved together with the following boundary conditions: (a) along the 
surface of the drop the air is saturated with vapor and hence n at the interface is the constant 
density of the saturated vapor n^, (b) far away from the drop the density approaches the 
constant ambient vapor density rzoo, and (c) the vapor cannot penetrate the substrate and 
hence dzU = at the substrate outside of the drop. Having found the vapor density, one can 
obtain the evaporation rate J = —D'Vn, where D is the diffusion constant. 

This boundary problem is mathematically equivalent to that of a charged conductor of 
the same geometry at constant potential if we identify n with the electrostatic potential and 
J with the electric field. Moreover, since there is no component of J normal to the substrate, 
we can further simplify the boundary problem by considering a conductor of the shape of 
our drop plus its reflection in the plane of the substrate in the full space instead of viewing 
only the semi-infinite space bounded by the substrate (Fig. |^. This reduces the number 
of boundary conditions to only two: (a) n = on the surface of the conductor, and (b) 
n = at infinity. The shape of the conductor (the drop and its reflection in the substrate) 
is now symmetric with respect to the plane of the substrate and resembles a symmetrical 
double-convex lens comprised of two spherical caps. This equivalent electrostatic problem of 
finding the electric field around the conductor at constant potential in the infinite space is 
much simpler than the original problem in the semi-infinite space. The reflection technique 
for finding the evaporation field on the basis of the analogy between the diffusion and the 
electrostatics was originally used by Deegan et al. j2Il 122] ■ 

Even in the circular geometry the equivalent problem is still quite complicated despite 
the visible simplicity. We consider an object whose symmetry does not match the symmetry 
of any simple orthogonal coordinate system of the three-dimensional space. In order to solve 
the Laplace equation, one is forced to introduce a special coordinate system (the so-called 
toroidal coordinates) with heavy use of the special functions. The full solution to this problem 
is provided in the Appendix. 

The evaporation rate depends only on the overall shape of the drop, and evaporation oc- 
curs in the same fashion from both phases. We assume that the evaporation is not influenced 
by any motion of the solute inside the drop, and the necessary amount of fluid can always 
be supplied to the regions of the highest evaporation near the contact line. Physically, high 
evaporation near the edge is what brings the solute to the contact line, and we assume that 
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Figure 5: Presence of particles in the deposit does not obstruct fluid evaporation at the 
edge of the drop. All the necessary fluid is supplied, and it is this motion of the fluid that 
brings the particles to the deposit phase. Also shown schematically is the fact that the 
boundary between the phases is vertical and the particles get stacked at full height between 
the substrate and the free surface of the drop. 



presence of the deposit does not obstruct the motion of the fluid (Fig. EI). Since the drop is 
thin and the contact angle is small, we will use expression 



for the evaporation rate (derived in the Appendix for no-solute drops in the limit <^ 1), 
which has the one-over-the-square-root divergence near the contact line intuitively expected 
from the electrostatics. The real situation may be different from the assumed above when p 
is large or comparable to 1, and the edge of the area where the evaporation occurs may be 
located near the boundary of the phases instead of the contact line. However, for small initial 
concentrations of the solute, the main order result will be insensitive to the exact location of 
the singularity of the evaporation rate: whether it is located at the contact line or near the 
boundary of the phases. We will further comment on this case of the "dry deposit" when we 
obtain the full system of equations. 

3 Principal equations 

Global conservation of fluid. The essence of the entire theory can best be summarized 
in one sentence: "It is all about the conservation of mass." Indeed, as we will see by the 
end of this section, all three governing equations obtained here represent the conservation of 
mass (or volume) in one form or another. 

We start from the global conservation of fluid in the drop. Since the amount of solute 
within the drop does not change during the drying process, the change of the entire drop 
volume is equal to the change of the amount of fluid. This fluid gets evaporated from the 
surface, and the total change of the fluid volume equals to the amount evaporated from the 
surface: 

dVL, = dV^ . (12) 

surf ^ ' 

By convention, superscripts F and S refer to the fluid and the solute components, respectively 
(while subscripts L and D continue to denote phases). The total change of the drop volume 



J(r) 



2 D(ns - rioo) 



(11) 
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is the sum of the volume changes of each phase: 



dV\,,, = dV, + dVn = ^d (^^^ , (13) 

where dVi and dV^ were found in the preceding section [Eqs. ^ and Q]- The volume of 
fluid evaporated from the surface can be determined from the known evaporation rate: 



dV^ 



dt 



surf 



^^'^ /rT(W2-rfr = -^^^!^^^:^, (14) 



P ' P 



where p is the fluid density. We neglected the gradient of h{r) with respect to unity (which 
is always legitimate for thin drops) and used J(r) of Eq. (fTT|) . Thus, Eqs. (fT^ . (fT^ . and 
(fT^ yield the first main differential equation of this section: 

^.d e\ IQDi 
dt \R J np 

This equation represents the global conservation of fluid in the drop and relates the time 
dependencies of 9(t) and R(t). 



Local conservation of mass. The next equation represents the local conservation of 
mass. There are two components in the liquid phase, and hence we write a separate equation 
for each of them. Since a free particle of the appropriate size reaches the speed of the flow 
in about 50 ns in water under normal conditions ^32', the solute particles are simply carried 
along by the flow, and the velocities of each component are identical at each point within the 
liquid phase [and equal to the depth-averaged fluid velocity v defined in Eq. ([T))]. The local 
conservation of fluid can be written in the form: 

V ■ [(1 - x)ih + i/)v] + ^yrT(W + dt [(1 -X){h + H)] = 0, (16) 

where x is the volume fraction of solute at a given point within the liquid phase, and each 
of the quantities {h + H), J, x, and v is a function of distance r and time t. [We drop the 
(V/i)^ part of the second term everywhere in this work since it is always small compared to 
unity for small contact angles.] This equation represents the fact that the rate of change of 
the fluid amount in a volume element (column) above an infinitesimal area on the substrate 
(third term) is equal to the negative of the sum of the net flux of fluid out of the column (first 
term) and the amount of fluid evaporated from the surface element on top of that column 
(second term); Fig. IHl illustrates the idea. A similar equation can also be written for the local 
conservation of solute, but without the evaporation term: 

V ■ [x{h + H)y] + dt [xih + H)] = 0. (17) 

Adding the two equations and employing the linearity of the differential operations, one 
obtains: 

V ■{{h + H)w] + - + dt{h + H) = {]. (18) 
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Figure 6: Conservation of mass: the liquid-vapor interface lowers exactly by the amount of 
fluid evaporated from the surface plus the difference between the outflow and the influx of 
fluid from the adjacent regions. 



This relation could have been obtained if we considered only one component with volume 
fraction 1 in the liquid phase, and this equivalence should be of no surprise: when the solute 
moves in exactly the same fashion as the fluid does, any differentiation between the two is 
completely lost (from the point of view of the conservation of volume). Note that if evapo- 
ration were too intensive, this equivalence would not hold, as there might be an insufficient 
amount of fluid coming into a volume element, and the solution could get completely dry 
(only the solute component would be left). We implicitly assume this is not the case for 
our liquid phase where the solute fraction is relatively small and the evaporation is not too 
strong. 

In circular geometry, due to the symmetry, the flow is radial and independent of (j). Thus, 
Eq. (jl8|) can be resolved with respect to the radial component of the velocity: 

^-(^' ^) = - ! m r i- + dth + dtH] rdr. (19) 
r[h + Hj Jo \p J 

Straightforward integration with h{r,t) of Eq. dH/dt of Eq. ©, and J(r) of Eq. (fTT|) 
and employment of Eq. (fT3j) for d{0 / R)/dt yield 



Vr{r,t) = 




This expression for the flow velocity at each point r within the liquid phase in terms of 
the time-dependent geometrical characteristics of the drop 9{t), R(t), and H(t) is a direct 
consequence of the local conservation of mass. 

With the velocity in hand, we can compute the time it takes an element of fluid initially 
located at distance from the center to reach the contact line. First, only the particles 
initially located near the contact line reach that contact line. As time goes by, the particles 
initially located further away from the contact line and in the inner parts of the drop reach 
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the contact line. Finally, the particles initially located in the innermost parts of the drop 
[i.e. near its center) reach the contact line as well. The more time elapsed, the more particles 
reached the contact line and the larger the area is where they were spread around initially. 
One can view this process as inward propagation of the inner boundary of the set of the 
initial locations of the particles that have reached the contact line by time t. As is easy 
to understand, the velocity of this front is equal to the negative of the vector of the fluid 
velocity at each point (the fluid and the particles move together towards the contact line 
while this front moves away from it, hence a minus sign). We label rj(t) the initial location 
of the solute particles that reach the phase boundary (and become part of the deposit ring) 
at time t. Since the solute particles from the outer areas of the drop reach the deposit phase 
sooner than the particles from the inner areas, this function is monotonically decreasing, and 
its derivative is simply related to Vr found in the preceding paragraph [Eq. (j^ ]: 

^ = -^rir,,t). (21) 
Thus, the second principal equation of this section is 

dvi _ 2D{ns - Uoo) Ri V'"'" ~ (^) 
dt vrp Tj w 

2 

This equation relates rj(t) to the time dependencies of the geometrical parameters of the 
drop [e{t), R{t), and H{t)]. 





Global conservation of solute. The volume of solute in the deposit phase V£ at time 
t is equal to the volume of solute outside the circle of radius rj(t) at time (since all the 
solute between rj(t) and Ri becomes part of the deposit by time t). The latter volume can 
be found by integrating h{r, 0) over the area swept by the fluid on its way from rj to the 
contact line and multiplying the result by the initial volume fraction of solute xf- 



Xi / h{r, 0) 2nrdr = V" 



1 - I — 

Rj 



(23) 



where V^"^ = nXiR^^^i/^ is the total volume of solute in the drop. On the other hand, the 
volume of solute in the deposit phase is just the constant fraction p of the volume of the 
entire deposit phase: 



V^^ =pVd. 



(24) 



Equating the right-hand sides of these two equations, taking the time derivatives of both 
sides, and making use of the already determined dVo of Eq. 0, we obtain the third principal 
equation of this section: 



XiRl^^i 

This equation represents the global conservation of solute in the drop. 



1- 4^ 



Ri 



d 
It 



ApHR—. 

dt 



(25) 
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Thus, we have four unknown functions of time: 6{t), R{t), H{t), and ri{t), and four 
independent differential equations for these functions: Eqs. (0), (fT3j) . (j^ . and (j^ . In 
reahty, we need only three of these functions: 9{t), R{t), and H{t); however, there is no 
simple way to eliminate rj(t) from the full system and reduce the number of equations. 
Having solved this system of equations, we will be able to fully characterize the dimensions 
of the deposit phase and describe the evolution of the deposit ring. The following section is 
devoted to the details and the results of this solution. 

Here we will only comment on how this system changes in the case of the completely dry 
solute p = 1. In this case there is no evaporation from the surface of the deposit phase, and 
the effective edge of the evaporating area is somewhere in the vicinity of the phase boundary. 
Assuming the same one-over-the-square-root divergence of the evaporation rate at r = i? 
instead oi r = Ri [which mathematically means substitution of R in place of Ri in Eq. (fTT|) ] 
and conducting a derivation along the lines of this section, one can obtain a very similar 
system of four differential equations. These equations would be different from Eqs. (0), (fT3j) . 
fl22|l ■ and ()25|) in only two minor details. First, Eqs. and ()22|1 would lose all indices i at all 
occasions of Ri {i.e. one should substitute R for all Ri in both equations). Second, p should 
be set to 1 in Eq. Apart from these details, the two systems would be identical. As we 
will see in the following section, this difference between the two systems is not important in 
the main order in a small parameter introduced below, and thus this "dry-solute" case does 
not require any special treatment contrary to the intuitive prudence. 

4 Results 

Analytical results in the limit of small initial concentrations of the solute. 

So far we have not introduced any small parameters other than the initial contact angle 
6i <^ 1. In particular, equations (0, (fT^ - (j^ . and were obtained without assuming 
any relation between p and Xi other than the non-restrictive condition Xi < P- order to 
find the analytical solution to this system, we will have to assume that Xi P- Then, we will 
solve the same system of differential equations numerically for an arbitrary relation between 
Xi and p. 

Assumption Xi ^P physically means that the solute concentration in the liquid phase is 
small — it is much smaller than the concentration of close packing or any other comparable 
number of the order of 1. This is the case for most practical realizations of the ring deposits 
in experiments and observations: the solute concentration rarely exceeds 10% of volume, and 
in most cases it is far lower. If the volume fraction of the solute is small, then the solute 
volume is also small compared to the volume of the entire drop. Hence, the deposit phase, 
which consists mostly of the solute, must also have small volume compared to the volume of 
the entire drop. Thus, if the initial volume fraction Xi is small, then the dimensions of the 
deposit ring must be small compared to the corresponding dimensions of the entire drop. 

Let us now introduce parameter e that is small when Xi/p is small. We do not fix its 
functional dependence on Xi/P for the moment: 




(26) 
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where / is an arbitrary increasing function of its argument. Then we postulate that the ring 
width is proportional to this parameter: 



R{t) =Ri[l- eW{t)\ , (27) 

where W{t) is an arbitrary dimensionless function and we explicitly introduced the dimen- 
sionality via Ri. Obviously, 1^(0) = 0. So far we simply wrote mathematically that the ring 
width is small whenever the initial volume fraction of the solute is small. Next, we introduce 
a dimensionless variable for the angle 0{t): 

Oit) = (28) 

where both 9{t) and 9i are small, while the newly introduced function 9{t) is arbitrary [in 
particular, ^(0) = 1]. Due to the geometrical constraint the height of the ring H{t) must 
be linear in small parameters e and 9i and directly proportional to the only dimensional scale 
Ri'. 

H{t) = e9iRiH{t), (29) 

where H{t) is yet another dimensionless function of time [H{0) = 0], related to functions W{t) 
and 9{t) by an expression similar to Eq. (jH]). The last dimensionless variable is introduced in 
place of the fourth unknown function rj(t): 

Vit) = l-{^J, (30) 

with the initial condition V{0) = 0. Thus, we introduced four new dimensionless variables in 
place of the four original ones and explicitly separated their dependence on small parameters 
e and 6'j. Finally, we define the dimensionless time r as: 

r=f, (31) 

where is a combination of system parameters with the dimensionality of time: 

tf = (32) 

^ 16D(n,-noo) ^ ' 

In the limit Xi/p ^ ^ this combination represents the time at which all the solute reaches 
the deposit phase; for finite Xi/p it does not have so simple interpretation. 

Substitution of all the definitions of the preceding paragraph into the original system of 
equations (0), (fT3j) . ((221), and and retention of only the leading and the first correctional 
terms in e yield the following simplified system of equations: 

d9 ~dW ~ d9 

-T- + eO-- - ?>eW— = -1, 34 
dr dr dr 
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dV 
dr 



AeW (V-^ 



2ev 



1 - eW (2\/-i - 1 



4eH 



dr 



Ae^H- 



^~dW 



l-eW 



(35) 



(36) 



p dr dr 

As is apparent from the last equation, parameter must be proportional to Xi/P- Since the 
separation of the ring width into e and W in Eq. (j^Tj) is absolutely arbitrary, parameter e is 
defined up to a constant multiplicative factor. Therefore, we set this factor in such a way 
that is equal to Xi/P- 



(37) 



This fixes the function / from the original definition ()26|) . 

The differential equations in the system (jnnj) - ()36|) are still coupled. However, in the main 
(zeroth) order in e, the equations clearly decouple: the second equation can be solved with 
respect to 9{t) independently of all the others, then the third equation can be solved with 
respect to V{t) independently of the first and the fourth, and finally the first and the fourth 
equations can be solved together as well. Thus, one can obtain the following main-order 
solution to the system of equations above with the appropriate initial conditions: 



fr) = 1 - r. 



V{r)= l-(l-r) 



3/4' 



2/3 



H{r) 
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(40) 
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dr'. 



Here B(a, b) 



X 



Jo 8H{t') (l-r')i/4 
dx is the complete beta-function, 62(0,6) = j^x^-'^il 



(41) 



xY~^ dx is the incomplete beta-function (a > 0, 6 > 0, and < ^ < 1), and the integral 
in the last equation cannot be expressed in terms of the standard elementary or special 
functions. In a similar fashion, systems of equations of the higher orders in e can be written 
[only the first-order corrections are kept in the system ()33|) -(|nU |) ]. and the higher-order terms 
can also be constructed up to an arbitrary order. 

A system of equations similar to our system was presented by Robert Deegan in 

works 121] . However, some terms of the first order in concentration were missing and no 
analytical solution to the system of equations was obtained in those works. Here we derive 
the equations in a systematic way. The analytical solution (| ^ -(PT |) is provided for the first 
time. 

How do our results ()38 p - ()41|) translate into the original variables? The first two of them 
[Eqs. (jHS|) and (jSHl)] reproduce earlier results. In terms of the dimensional variables Eq. (|HS|) 
represents the linear decrease of the angle between the liquid-air interface and the substrate 
at the phase boundary with time: 



e{t) = (1 - ^ 



(42) 
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Figure 7: Mass of a drying drop as a function of time. Experimental results, after Refs. [23 
1^ . The line running through the data is a linear fit. (Courtesy Robert Deegan.) 



[plotted by the solid line in Fig. IHt^a)]. This is a direct analog of Eq. (j64j) for the contact 
angle in the no-solute case, as is clear from the definition of [Eq. (j32p ]. So, angle 9 in the 
case of the finite-volume solute depends on time in exactly the same fashion as the contact 
angle in the no-solute case does. This expression also provides an interpretation of tf: it is 
the time at which the free surface of the liquid phase becomes fiat. Before tf this surface is 
convex, after it becomes concave and bows inward (until it touches the substrate). Thus, 
is generally not the total drying time. In the limit Xi/p ^ the height of the deposit ring 
is going to zero and the two times are the same. For finite values of this parameter the total 
drying time is longer than the time at which the liquid-air interface becomes fiat. Eq. (jl^ 
has been verified in the experiments [221 121] where the mass of the drop was measured as a 
function of time (Fig. [7|). Since the mass of a thin drop is directly proportional to 9, these 
results confirm the linearity of 9{t) during most of the drying process. 

The second equation (jH^ has a direct analog in the case of the zero-size solute particles. 
In the original variables it can be rewritten as 



which is identical to Eq. (3.24) of Ref. [3T] obtained for the zero-volume solute. Clearly, 
Ti = Ri when t = 0, and = when t = tf. According to Eqs. (j23|) and ([30|) . the fraction 




(43) 
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of solute in the deposit phase V^/V^ is 



Yl 



4/3 



(44) 



[plotted by the solid line in Fig. IHfb)]. This fraction is at t = and becomes 1 at t = tf. 
Thus, tf can also be interpreted as the time at which all the solute particles become part of 
the deposit phase. So far, the results of this finite-volume model coincide with the results of 
the zero- volume case considered earlier pTl I22j . 

However, the third and the fourth equations [Eqs. (j4(J|) - (j41|) ] represent entirely new results. 
In the dimensional variables they yield the height of the phase boundary H and the width 
of the deposit ring W = Ri — R, respectively: 



Hit) 



-OjRjH { — 



Wit) = Mr,w 



(45) 



(46) 



where functions Hir) and Wir) are plotted in Figs. |Hfc) and|Hl^d) (the solid curves). These 
results provide the sought dependence of the geometrical characteristics of the deposit ring 
on all the physical parameters of interest: on the initial geometry of the drop (i?i and 6i), 
on the initial solute concentration ixi), and on the time elapsed since the beginning of the 
drying process (t). If the time is considered as a parameter, they can also be used to obtain 
the geometrical profile of the deposit ii.e. the dependence of the height on the width), which 
we plot by the solid line in Fig. |H1 Note that the vertical scale of this plot is highly expanded 
compared to the horizontal scale since there is an extra factor of 6*^ <^ 1 in the expression for 
the height; in the actual scale the height is much smaller than it appears in Fig. El 

It is straightforward to obtain the asymptotics of Hir) and Wir) for early and late drying 
stages. At early times, both the height and the width scale with the drying time as a power 
law with exponent 2/3: 



[l + 0(r)] 




r 



«1) 



«1) 



(47) 



(4J 



P 

Thus, at early times H ^ OiW, which can also be deduced directly from Eq. (jHI) without 
obtaining the complete solution above. [The early-time exponent 2/3 was first obtained by 
Robert Deegan ^| without deriving the full time dependence (pn|) -PT |) .] At the end of 
the drying process, the height and the width approach finite values (which, apart from the 
dimensional scales, are universal, i.e. constants) and do so as power laws of (t/ — t) with two 
different exponents: 




Hil] 



(1-r) 
14^(1) 



7/4 



Oil - r)5/2 



(1 - r < 1) 



(49) 
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Figure 8: Results: dependence of the geometrical characteristics of the drop on time t. In 
each plot, the solid curve is the analytical result in the limit Xilv ~^ 0) while the other 
curves are the numerical results. Different numerical curves correspond to different initial 
concentrations of the solute; values of parameter Xilv ^-re shown at each curve, (a) Angle Q 
between the liquid-air interface and the substrate at the phase boundary, (b) Volume fraction 
of the solute in the deposit phase V^jV^ ■ (c) Height of the phase boundary H (in units of 
OiRi^ Xi/p)- (d) Width of the deposit ring W (in units of Ri\jxilv)- 
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Figure 9: Deposit ring profile: dependence of the height of the phase boundary H on the 
width of the deposit ring W . The solid curve is the analytical result in the limit Xilv ~^ 0; 
the other curves are the numerical results. Different numerical curves correspond to different 
initial concentrations of the solute; values of parameter Xilv ^'^'^ shown at each curve. The 
vertical scale is different from the horizontal scale by a factor of ^ 1. 
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where H{1) and W{1) are simply numbers: 




0.297, 



W{1) 
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(51) 



(52) 



Clearly, dH/dW = 0^(1 — r) and hence vanishes when r ^ 1. This fact can also be observed 
in the flattening of the analytical graph in Fig. IHl at late times. 

Dependence of the height and the width on the radius of the drop Ri, while intuitively 
obvious (since Ri is the only scale in this problem with the dimensionality of the length), has 
been verifled in experiments fl^ I24j . A linear flt has been obtained for the dependence of 
the ring width on the radius, in exact agreement with our flndings. 

Comparison to the experimental data for the dependence on the initial concentration 
of the solute is slightly less trivial. Our results predict that both the height H and the 

1/2 

width W scale with the initial concentration as Xi (at least, in the leading order for small 



concentrations). The same scaling prediction was also made by Robert Deegan 
However, his experimental results show a different exponent of xf- values 0.78 ± 0.10 and 
0.86 ± 0.10 were obtained for two different particles sizes (Fig. fTO|l . Why is the difference? 
The answer lies in the fact that the width measured in the experiments ^31 121] is not the 
full width of the ring at the end of the drying process, but rather the width of the ring at 
depinning. Depinning is a process of the detachment of the liquid phase from the deposit 
ring fFig. Illj). This detachment was observed experimentally in colloidal suspensions but has 
not been explained in full theoretically yet.^ An important observation, however, is that the 
depinning time {i.e. the time at which the detachment occurs and the ring stops growing) 
depends on the initial concentration of the solute. This dependence was also measured by 
Deegan (Fig. IT^ . and the resulting exponent was determined to be 0.26 ± 0.08. Thus, the 
width of the ring at depinning Wd scales with the initial concentration of the solute Xi 



Wdocx 




(xxf'Wlx 



,0.26±0.08 



(53) 



where td is the depinning time {td/tj oc Xi''^^^^"^)- As is apparent from Fig. ^1 the typical 
values of the depinning time are of the order of 0.4-0.8 tf. In this time range, function W{t) 
is virtually linear [the analytical curve in Fig. |H[d)]. Therefore, the dependence of Wd on Xi 
has the overall exponent of the order of 0.76 ± 0.08. It is now clear that both experimental 

^While the full explanation is yet to be developed, the naive reason for the depinning seems relatively 
straightforward. The pinning force depends only on the materials involved and is relatively insensitive 
to the value of the contact angle. At the same time, the depinning force is simply the surface tension, 
which is directed along the liquid-air interface and which increases as the contact angle decreases (since only 
the horizontal component of this force is important). Thus, the relatively constant pinning force cannot 
compensate for the increasing depinning force, and after the contact angle decreases past some threshold, the 
depinning force wins and causes the detachment. 
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Figure 10: Ring width normalized by the drop radius vs. initial concentration of the solute for 
two different particle sizes. Experimental results, after Refs. fl^ I24j . The two data sets are 
offset by a factor of 5 to avoid mixing of the data points related to the different particle sizes. 
The lines running through the data are linear fits in the double-logarithmic scale, which upon 
conversion to the linear scale yield power laws with exponents 0.78 ± 0.10 and 0.86 ± 0.10. 
(Courtesy Robert Deegan.) 



results 0.78 ± 0.10 and 0.86 ± 0.10 fall within the range of the experimental uncertainty of 
this approximate predicted value, and the theoretical dependence of the ring width on the 
initial concentration agrees with the experimental results quite well. 

Note that Robert Deegan j2Sl 12^ did not report direct measurements of the height of 
the deposit ring. The height was calculated from the data in hand, and thus the direct 
comparison to the experimental data is not available for the height. 

The square-root dependence of the height and the width on the concentration is in good 
agreement with general physical expectations. Indeed, the volume of the deposit ring is 
roughly proportional to the product of the height and the width. On the other hand, the 
height is of the same order of magnitude as the width since the ratio of the two is of the 
order of 6i (which is a constant). Thus, both the height and the width scale approximately 
as a square root of the ring volume. Finally, the volume of the deposit ring is proportional 
to the initial volume fraction of the solute: the more solute is present initially, the larger the 
volume of the deposit ring is at the end. Therefore, both the height and the width must scale 
as a square root of the initial volume fraction. It is rewarding that our complex calculation 
leads to the same results as this simple physical argument. 
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Figure 11: A photographic sequence demonstrating a depinning event. Experimental results, 
after Refs. 12^ ■ The view is from above, and the solid white band in the lower part of 
the frame is the ring; the rest of the drop is above the ring. The time between the first and 
the last frames is approximately 6 s; the major axis of the hole is approximately 150 /xm. 
(Courtesy Robert Deegan.) 
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Figure 12: Depinning time normalized by the extrapolated drying time vs. initial concentra- 
tion of the solute. Experimental results, after Refs. The line running through the 
data is a linear fit in the double-logarithmic scale, which upon conversion to the linear scale 
yields a power law with exponent 0.26 ± 0.08. (Courtesy Robert Deegan.) 
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Thus, the complete analytical solution to our model is available in the limit Xi/p ~^ O5 ^ind 
this solution compares favorably with the experimental data. Since the main-order solution 
in Xi/p is perfectly adequate, the difference between the original system of equations and 
the one for the "completely dry" case is not important. Indeed, the main-order results are 
identical in both cases, because one case in different from the other only by presence of R 
instead of Ri in a few places in the main equations, and this difference is of the correctional 
order in Xi/P- 

Numerical results for arbitrary initial concentrations of the solute. Apart from 
approaching the original system of equations (0), (fT3j) . (f^ . and analytically, we also 
solve it numerically. During this numerical procedure we do not presume that Xi/p is small, 
nor do we expand any quantities or equations in e or any other small parameters. Our main 
purpose is to reproduce the results of the first part of this section and to determine the range 
of validity of our analytical asymptotics. 

The typical values of Xi/p i^i most experimental realizations are of the order of 0.001-0.01, 
and thus, only the concentrations below approximately 0.1 are of practical interest. (Note that 
Xi/p = 0.1 corresponds to a quite substantial value of the small parameter e = a/cTT ~ 0.32.) 
Thus, we will concentrate on this range of Xi/p when describing the results despite the fact 
the numerical procedure can be (and have been) conducted for any ratio Xi/P- The general 
trend is illustrated well by the results in this range of concentrations. In the case of Xi 
comparable to p our model is not expected to produce any sensible results, as the entire 
separation of the drop into the two phases (the liquid phase and the deposit phase) is based 
on the assumption that the mobility of the solute is qualitatively different in the two regions. 
When Xi is comparable to p the two phases are physically indistinguishable, while the model 
still assumes they are different. 

We present our numerical results for the same quantities (and in the same order) as in 
our analytical results ()42|) and ()44|) - (jmj) . Since for arbitrary Xi/p time tf is not exactly the 
total drying time, there is a question of where (at what time) to terminate the numerical 
curves. By convention, we terminate all the curves in all the graphs at value of t/tf when 
all the solute reaches the deposit phase. In our model, it turns out that the time the last 
solute particles reach the deposit ring and the time the center-point of the liquid-air interface 
touches the substrate are about the same. For all the initial concentrations, the two times 
were numerically found to be within 0.1% of each other, and the curves are terminated at 
exactly this moment. Of course, in reality a small fraction of solute should stay in the 
liquid phase as long as the liquid phase exists, and so the moment the last solute particles 
reach the deposit phase should be after the moment the center-point touches the substrate; 
however, the amount of solute remaining in the liquid phase at touchdown is insignificant, 
and practically all the deposit has already formed. 

Numerical results for angle ^ as a function of time are shown in Fig.|Hfa). All curves behave 
almost linearly (as expected), however, the slope increases with concentration: formation of 
the ring in the drops with more solute finishes faster (in the relative scale of tf). The end 
of each curve demonstrates the value of the angle 6t at the moment the liquid-air interface 
touches the substrate. The analytical expression for this angle is dt = —2H/R for a thin drop. 
The absolute value of this angle increases with concentration, which is quite natural since for 
small concentrations the height of the ring grows as a square root of the concentration while 
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the radius of the hquid phase does not change substantially. Clearly, the numerical results 
converge to the analytical curve when x%lv ^ 0- 

Growth of the volume fraction of solute in the deposit phase V^jV^ with time is shown 
in Fig. lEfb) for various solute concentrations. This graph reconfirms the observation of the 
preceding paragraph that the solute transfer happens faster (in units of tf) for denser colloidal 
suspensions. All curves are terminated when volume fraction Vi)lV^ becomes equal to 1. 
(The apparent termination of the curve for Xilv = 0-1 earlier than that is an artifact of 
the plotting software.) As the corresponding analytical results do, the numerical plots of 
Figs. IH^a) and IHl^b) should presumably hold true independently of the geometrical details of 
the solute accumulation in the ring (which cannot be expected from the following plots for 
the ring height and width). 

The next two graphs represent the numerical results for the height [Fig.|Sl|^c)] and the width 
[Fig. IH^d)] of the deposit ring as functions of time. The ring profile, i.e. the dependence of 
the height on the width, is also shown in Fig.|Hl As the graphs depict, the ring becomes wider 
and lower (in the reduced variables) for higher initial concentrations of the solute. Since the 
volume of the ring is roughly proportional to the product of the height and the width, the 
decrease in height must be of the same magnitude as the increase in width. This can be 
qualitatively observed in the graphs. 

As a final piece of the numerical results, we create a double-logarithmic plot for the 
dependence of the height and the width on the initial concentration of the solute (Fig. IT!^ . 
The predicted square-root dependence on the initial concentration is seen to hold true for 
volume fractions up to approximately 10~^/^p for the height and up to approximately 10~^/^p 
for the width. The deviations for higher volume fractions are due to the increasing role of 
the correctional terms in e compared to the main-order terms represented by the solid lines. 
In this graph, as in all the results of this section, it is clear that our main-order analytical 
results provide an adequate description of all the functional dependencies in the range of the 
initial concentrations of experimental importance (0.001-0.01). 

In general, our numerical results complement and reinforce our analytical results, provid- 
ing a crosscheck of both methods. 

5 Discussion 

Both the analytical results of Eqs. ()42|) and ()44|) - ()46|) and the numerical graphs of Figs. |H1 and 
l^may be reproduced experimentally giving validation to the proposed model. Measurements 
of the profiles in Fig. El should be particularly easy to conduct (since there is no time depen- 
dence involved) and may confirm or refute the predicted robustness and the universality of 
the deposition profiles. 

While the main principles of the proposed model were laid down by Robert Deegan [23 
1^ . its analytical solution for small concentrations and its numerical solution for arbitrary 
concentrations are obtained here for the first time. The availability of the exact analytical 
solution demonstrated that the theoretical scaling of the deposit width at depinning with 
solute concentration indeed agrees with both measured values of the exponents. The earlier 
estimate of Refs. |^ 1^ based solely on the early-time exponent had an overlap with only 
one of the concentration exponents. Our numerical results also quantify the range of solute 
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Figure 13: Numerical results: log-log plot of the dependence of the height of the phase 
boundary H and the width of the deposit ring W on the initial volume fraction of the solute 
Xi- The main-order analytical results H oc \Jxi/v and W oc \J~Xifp are also provided for 
comparison. 



concentrations where the predicted square-root dependence of the width holds true. In gen- 
eral, Robert Deegan's results were not sufficient to obtain the proper scaling of the width 
with time anywhere beyond the very early drying stages; the results of this work provide 
that time scaling at all drying stages. All the presented results suggest that the deposit ring 
profile and its growth can be fully accounted for on the basis of the finite volume of the solute 
particles only and that the governing functional dependences are universal. 

One may notice that the curves in Fig. El end at some positive (non-zero) height. This 
indicates the solute is exhausted before the profile curves have a chance to return to the 
substrate, and the final shape of the deposit ring must have a vertical wall at its inner side. 
We believe this is an artifact of our model, which is inherently two-dimensional when fiows 
inside the drop are concerned. Thus, the vertical distribution of the solute was assumed 
homogeneous (the phase boundary is vertical and the particles get stacked uniformly at 
all heights), and we used the depth-averaged velocity (P) throughout this work. This is 
equivalent to assuming that vertical mixing is complete. This assumption is quite important, 
and the results are expected to get modified if the true three-dimensional velocity profile is 
used instead of the depth-averaged velocity. We expect that if a three-dimensional model 
were built and the dependence on z were taken into account for all the quantities then 
the discontinuous wall of the phase boundary would get smoothened and the height would 
continuously return to zero. A question remains whether such a model would be solvable 
analytically. 

Our model relies on the assumption that solute mobility is different in the so-called liquid 
and deposit phases. In essence, we assume that the mobility is in the deposit phase and 1 in 
the liquid phase. This assumption, while artificial in its nature, seems relatively reasonable 
when applied to this system. Indeed, in the physical situations near the close packing, the 
loss of mobility typically occurs over a quite narrow range of the concentration values, and 
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hence our assumption should work satisfactorily when the difference between Xi s^nd p is in 
the orders of magnitude. The higher the initial concentration is and the closer the two values 
are, the worse this assumption holds true and the more artificial the difference between the 
two phases is. Thus, the validity of any model based on this separation of the mobility scales 
decreases for higher initial concentrations of the solute. 

The model assumes that the free-surface slope between the liquid and the deposit phases 
is continuous. In fact, assumption (jH)) expressing this continuity is one of the four basic 
equations of this work. This assumption seems quite natural as well. Indeed, if the liquid 
is present on both sides of the phase boundary, the change in the slope of its free surface 
would cost extra energy from the extra curvature at the phase boundary, since the liquid-air 
interface possesses effective elasticity. Presence of this extra energy (or the extra pressure) at 
the location of the phase boundary is not justified by any physical reasons as all the processes 
are slow and the surface is in equilibrium. In equilibrium, the surface shape must have 
constant curvature past the phase boundary since the entire separation into the two phases 
is quite artificial as discussed above. Presence of the particles below the liquid-air interface 
does not influence the surface tension, and thus the liquid surface (and its slope) should be 
continuous at the phase boundary. If the density of the particles matches the density of the 
liquid (which was the case in the experiments), nothing prevents the particles from filling up 
the entire space between the substrate and the liquid-air interface, thus providing the growth 
of the upper edge of the deposit phase along the liquid-air interface. This is particularly true 
for the thin drops discussed here, where vertical mixing is intensive, where the free surface 
is nearly horizontal, and where the problem is essentially two-dimensional. However, the 
equality of the slopes on both sides of the phase boundary does not seem inevitable, and 
one may think of the situations when it does get violated. One example might be the late 
drying times, when the deposit growth is very fast [Fig. Efb)] and hence the deposition may 
occur in some non-regular manner inconsistent with this slow-process description. Other 
examples may be related to gravity (slightly unequal densities of the particles and the fluid) 
or convection. This assumption can possibly be checked experimentally, and if condition (jH)) 
is found violated, an equivalent constraint dependent on the details of the deposit-growth 
mechanism must be constructed in place of Eq. (jH} . 

Another inherent assumption of our model is related to the evaporation rate J{r). Pres- 
ence of the solute inside the drop was assumed not to affect the evaporation from its surface. 
This is generally true when the evaporation is not too fast and the deposit phase is not too 
thick and not too concentrated. When these conditions are not obeyed, presence of a thick 
or concentrated layer of the solute on the way of the liquid moving from the phase boundary 
to the contact line may create a strong viscous force. This viscous force would prevent the 
necessary amount of fluid from being supplied to the intensive-evaporation region near the 
contact line. Generally, we assumed throughout this work that the viscous stresses are not 
important, and this is valid whenever v <^ a/Sr]. In the deposit phase, the velocity is large 
due to the proximity to the contact-line divergence of the evaporation rate, and the effective 
viscosity is large due to the high concentration of the solute. Thus, this condition may get 
violated and the viscosity may become important in the deposit phase, slowing down the 
supply of the liquid and ultimately making the deposit dry. Obviously, this affects the evap- 
oration rate, and the functional form of the evaporation profile changes. Simple assumption 
that the evaporation rate stays of the same functional form, but with the divergence at the 



27 



phase boundary (at R) instead of the contact hne (at i?,), was shown above not to affect our 
main-order results. Thus, our resuhs appear to be relatively insensitive to the exact loca- 
tion of this divergence within the (narrow) deposit phase. (In reality the evaporation edge 
would be somewhere between the original contact line and the phase boundary, i.e. the real 
situation is intermediate between the two considered.) However, the deposit could modify 
the evaporation rate J(r) in other ways. When there is a dry deposit ring just outside the 
liquid phase, the entire functional form of J may change, and the Laplace equation for an 
equivalent electrostatic problem must be solved anew with additional boundary conditions 
responsible for the presence of the dry solute rim and the modified evaporation at the edge. 
As we already saw in the Appendix, this is the most complicated part of the problem, and 
the mathematics can become prohibitively complex. Thus, finding the exact form of J may 
be a formidable task. One way around is in creating such evaporating conditions that the 
functional profile is simpler, for instance, J is just a constant. This would be more difficult 
to control experimentally, but would be much easier to treat analytically. The unavailability 
of the exact analytical form for J seems to be the biggest open question in this class of 
problems [23 EH- 

The equilibrium surface shape of the liquid phase is a spherical cap (0). This is a rigorous 
result valid during most of the drying process. However, when h{0, t) becomes negative 
and exceeds H(t) in its absolute value {i.e. when the center-point touches the substrate), the 
surface shape is no longer spherical. Moreover, a new element of the contact line is introduced 
in the center of the drop in addition to the original contact line at the perimeter, and the entire 
evaporation profile gets modified in addition to the modified surface shape, thus infiuencing 
all the other quantities. Our treatment does not account for the small fraction of the drying 
process occurring after this touchdown (which is a change in topology of the free surface, 
and thus requires a separate treatment after it happened). First of all, the amount of liquid 
remaining in the drop at this moment is of the order of e compared to the original volume, 
and hence it would not modify our main-order analytical results. Second, as our numerical 
calculations show, at touchdown practically all the solute is already in the deposit phase, and 
the remaining amount of solute in the liquid phase is insignificant. Thus, within our model, 
the remainder of the drying process cannot modify the deposit ring substantially, and hence 
this neglect of the late-time regime seems well justified. Experimentally, the inner part of the 
deposit ring is different from our prediction (which is a vertical wall) and appears to have a 
spread shelf. Presence of this tail in the deposit distribution can be caused by several features 
absent in our model. Its inherent two- dimensionality may be one of these shortcomings (as 
discussed above); the account for the dynamical processes occurring after the deposit phase 
has already been formed {e.g. avalanches of the inner wall) may be another missing feature. 
Absence of the treatment of the late-time regime may be among these reasons infiuencing 
the final distribution of the deposit as well. A more detailed account for the effects of this 
late-time regime might be required in the future. 

This work was completed as a part of the Ph.D. dissertation research supervised by Thomas 
A. Witten. The author acknowledges valuable input from Todd F. Dupont and Robert R. Deegan. 
This work was supported in part by the National Science Foundation MRSEC Program under award 
number DMR-0213745. 
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Appendix: On the evaporation rate 



The purpose of this section is to obtain the evaporation rate from the free surface of a round 
sessile drop on the substrate. Since presence of the solute is irrelevant to this purpose (at 
least for its low concentrations), one may assume that the solute is simply absent and the 
drop is just pure water. We first consider the generic problem with an arbitrary contact angle 
6', and then find the appropriate limit of interest ^ <^ 1. If the radius of the drop footprint 
on the substrate is Ri, then its surface shape for an arbitrary contact angle is given by 



hir,t) 



- R,coie{t). (54) 



Despite the fact this problem is two centuries old, some results are presented here in their 
closed analytical form for the first time, and some correct earlier expressions. 

Our task involves solution of the equivalent electrostatic problem (the Laplace equation) 
for the conductor of the shape of the drop plus its reflection in the plane of the substrate (kept 
at constant potential, as a boundary condition). In the case of the round drop the shape of 
this conductor resembles a symmetrical double-convex lens comprised of two spherical caps. 
The system of orthogonal coordinates that matches the symmetry of this object (so that 
one of the coordinate surfaces coincides with the surface of the lens) is called the toroidal 
coordinates (a,/?, 0), where coordinates a and j3 are related to the cylindrical coordinates r 
and z by 

Ri sinh a Ri sin (3 

cosh a — cos (3 ' cosh a — cos (3 ' 

and the azimuthal angle has the same meaning as in the cylindrical coordinates. Solution to 
the Laplace equation in the toroidal coordinates involves the Legendre functions of fractional 
degree and was derived in a book by Lebedev jSHj. The electrostatic potential or vapor 
density is independent of the azimuthal angle (f) and reads 

n(a, j3) = Tioo + ijis — noo) \/2(cosh a — cos j3) x 



X / _^ ^ ^_ — ;^^P_i/2+iT-(cosha) dr. (56) 



cosh 9t cosh(27r — (3)t 
cosh vrr cosh('7r — 9)t 

Here rig is the density of the saturated vapor just above the liquid-air interface (or the 
potential of the conductor), rioo is the ambient vapor density (or the value of the potential at 
infinity), and are the Legendre functions of the first kind (they are real valued). 

The surface of the lens is described by the two coordinate surfaces (3i = n — 9 and P2 = 7^ + 9, 
and the /3 derivative is normal to the surface. The evaporation rate from the surface of the 
drop is therefore given by 

J (a) = df3n{a, /3) l^^^n+p^ = D ^^^^^^ dpn{a, p) l^^g^.g , (57) 

where D is the diffusion constant and h/^ = i?j/(cosha — cos/5) is the metric coefficient 
in coordinate f3. (Note that an incorrect expression for J with a plus sign in the metric 
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coefficient was used in Eq. (A2) of Ref. |2I|-) Thus, the exact analytical expression for the 
absolute value of the evaporation rate as a function of r is available: 



D{ns - rioo) 

J{r) = B 

Ki 

r°° cosher 



- sin e + V2 (cosh a + cos Of^^ x 
2 



rco cosh uT 

/ tanh [(tt — 9)t] P_i /2+ir (cosh a) rdr 

Jo cosh TTT 



(58) 



where the toroidal coordinate a is uniquely related to the polar coordinate r on the surface 
of the drop: 

Ri sinh a ^^^^ 
cosh a + cos 9 

Expression ()58|) is valid for an arbitrary contact angle 6 and corrects an earlier expression 
of Ref. [2Z] [Eq. (28)] where a factor of in the second term inside the square bracket is 
missing. 

The expression for the evaporation rate is not operable analytically in most cases, as it 
represents an integral of a non-trivial special function (which, in its turn, is an integral of some 
simpler elementary functions). In most cases, it is necessary to recourse to the asymptotic 
expansions in the contact angle 6 in order to obtain any meaningful analytical expressions. 
However, there is one exception to this general statement. An important quantity is the 
total rate of water mass loss by evaporation dM/dt, which sets the time scale for all the 
processes. This total rate can be expressed as an integral of the evaporation rate (defined as 
the evaporative mass loss per unit surface area per unit time) over the surface of the drop: 

dM r _ / — — , . , fR 



J(r)^l + {Vhfrdrdcp = - j J{r)^fl+~{d^ 27irdr, (60) 



dt J A 

where the ffist integration is over the substrate area A occupied by the drop. This expression 
actually involves triple integration: one in the expression above as an integral of J(r), another 
in the expression for J(r) as an integral of the Legendre function of the ffist kind, and the third 
as an integral representation of the Legendre function in terms of the elementary functions. 
However, it is possible to simplify the above expression significantly and reduce the number 
of integrations from three to one. Investing some technical effort and using Eq. (2.17.1.10) 
of Ref. one can obtain a substantially simpler result that does not involve any special 
functions at all: 



—— = -nRiD{ns - 
dt 



sin 9 



1 + cosh26'r , , 

+4 / ; tanh vr - 9)t] dr 

sinh27rr ' ^ 



I + COS 9 

^ (61) 



(not reported in the literature previously). This result together with the expression for the 
total mass of water 

/■^',/ , „3Cos3^-3cos^ + 2 

M = p h{r, t) 27rrdr = npR^ — ^- 62 

JO 3 sm 9 

(where p is the water density) provides a direct method for finding the time dependence 
of 9 for an arbitrary value of the contact angle. Combining the time derivative of the last 
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expression with result one can obtain a single differential equation for 6' as a function 
of time t: 

de D{ns - n^) ^ T sin6' 

— (l + cosb') 



dt pR: 



2 



1 + cos 61 



r°°l + cosh26'r , , , 

+4 / tanh (tt - 6')r dr 

Jo sinh 2ttt 



(63) 



Having determined the dependence 6{t) from this equation, one can obtain the time depen- 
dence of any other quantity dependent on the contact angle, for instance, the time dependence 
of the mass from relation (jU^ . or any other geometrical quantity considered earlier. 

In practice, however, the analytical calculations in a closed form cannot be conducted 
any further for arbitrary contact angles, and we will use the limit of small contact angles in 
all the subsequent analytical calculations. Besides being the limit of our interest and most 
practical importance, this hmit is also perfectly adequate even for quite substantial angles, 
as will be seen in a moment. 

Expanding the right-hand side of Eq. ()63p in small 6, we immediately obtain that the 
contact angle decreases linearly with time in the main order of this expansion: 

e = e,(^i-^y (64) 

where we introduced the total drying time tf defined in terms of the initial contact angle 

tf = (65) 

In the main order, the total rate of water mass loss is constant and the water mass also 
decreases with time linearly: 

This linear time dependence during the vast majority of the drying process was directly 
confirmed in the experiments [23123; see Fig.[7| The dependence of the evaporation rate (pT|) 
on radius (linearity in Ri) was also confirmed experimentally and is known to hold true for 
the case of the diffusion-limited evaporation pij. 

In Fig. HH we plot the exact numerical solution for M{t) based on Eqs. (jHSl) and 
for several values of the initial contact angle 9i together with the small-angle asymptotic of 
Eq. In this figure. Mi is the initial mass of water in the drop defined by the pre-factor 

in Eq. (jMj) . [Note that t/ is not the total drying time for each 9i\ instead, it is just the 
combination of the problem parameters defined in Eq. (j65|) . which coincides with the total 
drying time only when 9i 0.] Fig. El demonstrates that the small-angle approximation 
works amazingly well up to the angles as large as 45 degrees, and therefore, no precision or 
generality is lost by working in the limit of small contact angles for the typical experimental 
values of 6',. Lastly, we note that the large- angle corrections may be responsible for the 
observed non-linearity of the experimentally measured dependence M(t), as is clear from the 
comparison of Fig. El (theory) and Fig. [3 (experiment). 
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Figure 14: Numerical results: dependence of water mass M on time t. Different curves 
correspond to different initial contact angles; values of parameter 9i are shown at each curve. 
The analytical result [Eq. (jMj) ] in limit ^ is also provided (the solid curve). 



Expression for the evaporation rate ()58|1 becomes particularly simple in the limit of small 
contact angles. Employing one of the integral representations of the Legendre function in 
terms of the elementary functions (Eq. (7.4.7) of Ref. |25j), it is relatively straightforward to 
obtain the following result: 

Ji.n = B cosh- (^^0), (67) 

which, upon identification cosh a = {R^ + r'^)/{Rf — r^) for 6 = 0, can be further reduced to 
Eq. (jllj) . Thus, for thin drops the expression for the evaporation rate reduces to an extremely 
simple result featuring the one-over-the-square-root divergence near the edge of the drop. The 
same result could have been obtained directly if we solved an equivalent electrostatic problem 
for an infinitely thin disk instead of the double-convex lens. It is particularly rewarding that 
after all the laborious calculations the asymptotic of our result is in exact agreement with 
the predictions of a textbook (see Ref. jHSj for the derivation of the one-over-the-square-root 
divergence of the electric field near the edge of a conducting plane in the three-dimensional 
space). Eq. (fTT|) is the result we were looking for in our case of the thin circular drops. 

For the sake of completeness, it is also interesting to note the opposite limit of the ex- 
pression ()58|) . when the surface of the drop is a hemisphere (6' = 7r/2). In this limit, a similar 
calculation can be conducted, and the uniform evaporation rate is recovered: 

J{r) = ^^^^^p^ i9^n/2). (68) 

This result is also in perfect agreement with the expectations; the same result could have 
been obtained if we directly solved the Laplace equation for a sphere (the hemispherical 
drop and its reflection in the substrate). The uniform evaporation rate is a result of the full 
spherical symmetry of such a system. Similar exact results can also be obtained for a few 
other discrete values of the contact angle {e.g. for 9 = 7r/4). 
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